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1. INTRODUCTION 

It has long been recognized 1 that spray characteristics and droplet vaporization physics play a fundamental role in 
determining the stability behavior of liquid fuel rocket motors. Such instabilities are generated by in-phase 
interactions of rocket acoustics with combustion fields. This paper gives some results of numerical parametric 
studies conducted to investigate the interactions of reacting and non-reacting sprays with acoustic waves. The results 
were obtained using a multi-dimensional algorithm which couples gas and liquid phase physics. The algorithm 
includes a new treatment for spray injection through direct solution of equations for the spray producing liquid core. 
A brief outline of the methodology is given in the next section. 

2 . TECHNICAL DISCUSSION 

2.1 General Approach 

The mathematical formulation provides for solution of the interaction between the liquid spray droplets and the 
surrounding gas. The formulation used in this study is based on a Lagrangian/Eulexian approach. Here, gas phase 
behavior is treated on a fixed (Eulerian) mesh. Behavior of liquid droplets such as heating, evaporation, and 
displacement is treated in a Lagrangjan manner. Interaction physics is contained in the source terms evaluated during 
the droplet-phase calculations, and subsequently introduced into the gas-phase formulation. 

Droplets are introduced into the computational domain following a separate computation of the liquid jet Liquid jet 
equations are solved on a separate grid embedded in the main computational domain. The jet motion and shape 
depends on a separately specified mass stripping rate, a semi-empirical formulation simulating the atomization 
process. Droplet sizes are specified using a semi -empirical expression based on linear stability theory. In order to 
conserve computational effort, droplets are treated as a “droplet cloud” with a number density dependent on the droplet 
mass and the stripping rate. The liquid droplets are allowed to evaporate and provide a source of fuel to the 
combustion calculations. 

The final results for liquid and gas phases are obtained by combining the Eulerian and Lagrangian calculations using 
the fractional step technique with operator splitting. Details of the model are given below. 

Gas Dynamics 

Gas phase equations of motion can be written in the general conservation form: 

$ + V • F c - “V * Fj = ^ + Sow.* ^ 

where U , F c , and Fj represent the conserved variables, and the convective and diffusive fluxes, respectively. The 
c * 

source term vectors “P and 0 ""* are due to droplet evaporation and combustion processes, respectively. 

Equation (1) is solved using the method of fractional steps where each of the convective, diffusive and source fluxes 
are solved for separately. The fractional steps are. 
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If each fractional step is solved using a second order scheme, second order accuracy erf the complete algorithm can be 
maintained by using Strang-type operator splitting. 2 Thus if Lg, Lq, L^, and Lcomb are the individual solution 
operators which advance the solutions through a time increment At, the solution at new time level can be expressed 
as. 


ir ' 2 = (L?L%LlL'l%LLi 
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(3) 


Such a solution procedure is widely used in combustion simulations; examples of its application can be found in 
Oran and Boris. 3 Note that first order solution requires only one sweep of the solution operators. 

2.2 Liquid Phase Equations 


Solution of liquid phase equations forms one fractional step on the overall computational scheme. These equations 
themselves are solved in fractional steps; the two main ones being formation and transport of liquid jets and spray 
droplets. 


Liquid Core 

The liquid core equations are solved on a separate mesh embedded in the main computational domain. Assuming an 
incompressible flow of liquid with negligible variations in velocity across the thickness, the equations describing the 
motion of an isothermal liquid core can be cast in the standard conservative form given by: 


au . . 

~3T + ~dx~ S ' 



(4) 


Note that the v-momentum equation is independent of the solution fat h or u. 

Hie source vector $ contains terms due to mass loss ^ caused by stripping of droplets, as well as external forcing 
functions / and g. Function / represents the effects of axial viscous drag, and g accounts for acoustic forcing caused 
by sound waves impinging on the liquid/gas jet interface. This formulation temporarily ignores the effects of surface 
tension and the axial pressure gradient, and introduces a mass stripping rate * to account for losses due to droplet 
production. Both the stripping rate and the resultant droplet diameters are specified using semi-empirical correlations 
which depend on the local relative velocity, gas and liquid densities, and surface tension. 4 Figure 1 below illustrates 
the development of a typical jet 



Figure 1. Development Time-History of a Typical Jet The succeeding profiles are separated by SO times teps. 



Since Equation (4) is hyperbolic, it may be solved using techniques similar to those employed in the gas-phase 
fractional step. However, since by ignoring the axial pressure gradient only one characteristic exists for the system, 
the solution is much simpler. 

Given droplets produced by the liquid core, droplet physics are computed using the spray modeL 

Spray Model 

The spray model used in this study is based on the Particle-In-Cell (PIC) formulation first proposed by Crowe et at 5 
and used with some variations by Moldavi and Sirignano. 6 In this approach, the droplet motion and evaporation are 
found by integrating the individual droplet Lagrangian equations of motion. Determination of source terms denoting 
droplet-gas exchange of mass, momentum, and energy is subsequently obtained from general balance equations at 
each computational cell. An extended form of this model accounts for 

• viscous and pressure drag on a droplet, including droplet deformation; 

• quasi-steady droplet heatup; and 

• quasi-steady vaporization. 


This fractional step convects and vaporizes the individual droplets and calculates the gas/liquid interaction source 
terms which are used in gas-phase calculations. This completes the description of the modeL 

3 . RESULTS 


Numerical experiments were performed to determine the effects of various parameters on spray characteristics and 
combustion response. The parameters included droplet size, stripping rate, spray angle, and chemical kinetics. A 
sample of the findings is given below. 


3.1 Effect of Droplet Site 


Droplets are inserted into the computational domain at the left, with initial velocity of 10 m/s. Transverse velocity 
is chosen at random to give a specified mean spray angle. Figure 2 shows unreacting spray shape and penetration for 
two drop at time t = 0.02s. This corresponds to droplet transit time in an in viscid flow. Note the difference in 
penetration length of the two sprays. __ 



(a) 


<b) 


Figure 2. Penetration into Quiescent Medium of Spray Composed of (a) 150p drops; and 9b) 1500ji drops 
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Liquid core is introduced into the domain at left The core develops in time until it reaches steady state. At each 
time step, droplets are stripped from the core and introduced into the domain. Initial drop velocity is that of the core 
with a random transverse component selected for a specified spray angle. Spray shapes are given in Figure 3 for a 
100 m/s jet producing 70|t drops, at time t = .001 s, for three spray angles. 



(a) (b) (c) 


Figure 3. Spray Shape and Penetration for Spray with a Liquid Core for Spray Angles of (a) 2°; (b) 4°; and (c) 6° 
Note the presence of the characteristic bulge. In the above calculations, the core extended to .06 m. 


This simulation is similar to that given in 3.2 above, except that the domain is now treated as a resonant tube with a 
20% cosine initial pressure perturbation, and the spray is allowed to react. Two reaction mechanisms for hydrogen- 
oxygen reaction: (a) instantaneous kinetics; and (b) (me step finite rate kinetics. In this simulation, liquid oxygen is 
injected into a hydrogen atmosphere. During the injection, the initial pressure inhomogeneity in the gas phase 
repeatedly traverses the domain, forms a shock and interacts with the spray. Figures 4 shows the effect of shock 
interaction of spray and the liquid core. 



Figure 4. Spray and Associated Liquid Core Profiles for (a) non-reacting spray; (b) reacting with infinite rate 

chemistry; and (c) reacting with one step finite rate model 
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Note the large spray deformations in the reacting case. These deformations are caused by the compressive effects of 
shock impinging on the spray coupled to the combustion generated expansion process. The figure also indicates the 
presence of disturbances on the surface of the liquid core. 

Figure 5 shows some typical end-wall pressure traces. Along with Figure 4, this figure shows that there is little 
difference in results for the two simple kinetic models used in this study. A more sophisticated chemistry model is 
needed to quantify the combustor sensitivity to chemical kinetics. 




(d) 


Figure 5. End Wall Pressure Time Histories: (a) base flow (no spray); (b) cold spray; (c) instantaneous 

kinetics; and (d) one step finite rate kinetics 


4. CONCLUSIONS 

A coupled numerical algorithm has been developed for studies of c omb ustion instabilities in spray -driven liquid 
rocket engines. The model couples gas and liquid phase physics using the method of fractional steps. Also 
introduced is a novel, efficient, methodology for accounting for spray formation through direct solution of liquid 
phase equations. Pre liminar y parametric studies show marked sensitivity of spray penetration and geometry to 
droplet diameter, considerations of liquid core, and acoustic interactions. Less sensitivity was shown to the 
c ombustion model type although more rigorous (multi-step) formulations may be needed for differences to become 
apparent Numerical studies are now planned to fully quantify these effects. 
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